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Abstract In this work, we calculate the self-similar longitudinal velocity correlation function and 
the statistical properties of velocity difference using the results of the Lyapunov analysis of the fully 
developed isotropic homogeneous turbulence just presented by the author in a previous work pQ. There, 
a closure of the von Karman-Howarth equation is proposed and the statistics of velocity difference is 
determined through a specific analysis of the Fourier-transformed Navier-Stokes equations. 

The correlation functions correspond to steady-state solutions of the von Karman-Howarth equation 
under the self-similarity hypothesis introduced by von Karman. These solutions are numerically deter- 
mined with the statistics of velocity difference. The obtained results adequately describe the several 
properties of the fully developed isotropic turbulence. 

Keywords Self-Similarity • Lyapunov Analysis • von Karman-Howarth equation • Velocity difference 
statistics 

PACS 47.27.-i 



1 Introduction 

A recent work of the author dealing with the homogeneous isotropic turbulence [T], suggests a novel 
method to analyze the fully developed turbulence through a specific Lyapunov analysis of the relative 
motion of two fluid particles. The analysis expresses the velocity fluctuation as the combined effect of 
the exponential growth rate of the fluid velocity in the Lyapunov basis, and of the rotation of the same 
basis with respect to the fixed frame of reference. The results of this analysis lead to the closure of the 
von Karman-Howarth equation and give an explanation of the mechanism of the energy cascade. A 
constant skewness of the velocity derivative du r /dr is calculated which is in agreement with the various 
source of the literature. Moreover, the statistics of the velocity difference can be inferred looking at the 
Fourier series of the velocity. This is a non-Gaussian statistics, where the constancy of the skewness 
of du r /dr implies that the other higher absolute moments increase with the Taylor-scale Reynolds 
number. 

The present work represents a further contribution of Ref . [if] . Here, the self-similar solutions of the 
von Karman-Howarth equation are numerically calculated using the closure obtained in the previous 
work and the statistics of the velocity difference is obtained. 
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2 Analysis 

For sake of convenience, this section reports the main results of the Lyapunov analysis obtained in the 
Ref. pp. 

As well known, the pair correlation function / of the longitudinal velocity u r for the fully developed 
isotropic and homogeneous turbulence, satisfies the von Karman-Howarth equation [2J 



dt u 2 \ dr 2 r dr J dr 2 

the boundary conditions of which are 

/(0) = i,aM = o 

* (2) 
lim /(r) = 

r — >oo 

Into Eq. ((T|), r is the separation distance and u is the standard deviation of u r , which satisfies [2J, [3] 
du 2 o d 2 f , . 

This equation gives the rate of the kinetic energy and is determined putting r = in the von Karman- 
Howarth equation [2], [3]. The function K(r), related to the triple velocity correlation function, repre- 
sents the effect of the inertia forces and expresses the mechanism of energy cascade. In accordance to 
the Lyapunov analysis presented in Ref [I], the analytical expression of K (r) is in terms of / and on 
its space derivative 



The proposed expression of K{r) satisfies the two conditions dK/dr(0) = 0, K(Q) = 0, which represent, 
respectively the condition of homogeneity and the fact that K does not modify the fluid kinetic energy 
[2] , [3] . The skewness of the longitudinal velocity difference Au r is calculated as [3] 

Hl(T) = _ m (5) 

((Au r r-f 2 (2(1 - /(,-))) 3/2 
where k{r) is the longitudinal triple velocity correlation function, related to K{r) through 

K[r) - u 3 + -) k(r) (6) 



dr 



r 



The result of this analysis is that the skewness of du r /dr is a constant which does not depend on 
the Reynolds number, whose value is i?3(0) = —3/7 pQ. The solutions of the von Karman-Howarth 
equation provide second and third dimensionless statistical moments of Au r . In line with the analysis 
of Ref. PQ, the higher moments are consequentely determined, taking into account that the analytical 
structure of Au r is expressed as 

Au r _ £ + ^(x(^ 2 -l)-(C 2 -l)) 



Equation (J7J) , which arises from statistical considerations about the Fourier-transformed Navier-Stokes 
equations in fully developed turbulence, expresses the internal structure of the isotropic turbulence, 
where £, rj and £ are independent centered random variables which exhibit the gaussian distribution 
functions and p(C) whose standard deviation is equal to the unity. Thus, the moments of 

Au r are easily calculated from Eq. ([7]) lj 

s iSy^ -jTTwh^U*) " kl( " w - " - ie - 1))l> (8) 
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where 

k 

((x(^ 2 -i)-(c 2 -i)) fe ) = E 



2 = 



(-xr<(c 2 -i)w -ir-*> 



(9) 



(=0 



In particular, the third moment or skewness, H3, which is related for the energy cascade, is 
H = 8^ 3 (x 3 - 1) 
3 (l + 2V> 2 (i + X 2 )) 3/2 

As far as ip is concerned, this is a function of the separation distance and of the Reynolds number pQ 



(10) 



ip(r,R) = 



R 



15V15 



(11) 



where i? = uX^jv and Ay = u/ ' \J {(du r / dr) 2 ) are, respectively, the Taylor-scale Reynolds number and 
the Taylor scale, whereas the function ip(r) is determined through Eq. (JTUJ) as soon as H 3 (r) is known. 
The parameter \ is a l so a function of R which is implicitly calculated putting r = into Eq. (|10|). i.e. 

m 



8V>o 3 



(1-x 3 ) 



(1 + 2VV(1 + X 2 )) 3/2 



3 
7 



(12) 



with = VK-R, 0) and ip = 1-075 [T]. 

From Eqs. Q and (fTT|) . all the absolute values of the dimensionless moments of Au r of order greater 
than 3 rise with i?, indicating that the intermittency increases with the Reynolds number. 
The PDF of Au r can be formally expressed through distribution functions p(£), p(i]) an d p(C)j using 
the Frobenious-Perron equation 



F(Au' r )= / p(£)p(rj)p(()5(Aur-Au' r )dZdrid{ 

J£ Jrj JQ 

where 5 is the Dirac delta. 

The spectrums E(k) and T(n) are the Fourier Transforms of / and K [3], respectively 



(13) 



'E(k)' 


1 r°° 


~u 2 f{r)~ 


T{n)_ 


T Jo 


K{r) 



K 2 r 2 



sin nr 



cos Kr I dr 



(14) 



3 Self-Similarity 

An ordinary differential equation for describing the spatial evolution of / is derived from the von 
Karman-Howarth equation under the hypothesis of self-similarity and using the closure given by Eq. 
((4]). This equation represents a boundary problem which is then transformed into an initial condition 
problem in the variable r. 

Far from the initial condition, it is reasonable that the simultaneous effects of the mechanism of the 
cascade of energy and of the viscosity act keeping / and E(n) similar in the time. This is the idea of 
self-preserving correlation function and turbulence spectrum which was originally introduced by von 
Karman (see ref. [4j and reference therein). 

In order to analyse this self-similarity, it is convenient to express / in terms of the dimensionless 
variables f — r/Xr and t = tu/\r, i-e., / = f(t, f). As a result, Eq. (fTJ) reads as follows 



df_ Xrd_ 

di u dt 



tu 
At 



df f dXr 
dr u dt 



1-fd.f 



2 df R \ df 2 f df 



R' 



(15) 
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where d 2 f/dr 2 (0) = — 1. This is a non-linear partial differential equation whose coefficients vary in 
time according to the rate of kinetic energy 

du 2 10vu 2 . . 

-& = -^r (16) 

If the self-similarity is assumed, all the coefficients of Eq. (fl"5l) must not vary with the time [2], [3], 
thus one obtains 

R = const (17) 

At d ( tu\ , . 

a\ = - — = const (18 

u dt \ At J 



1 gJAt 

ai = — = const, (19) 

u dt 

From Eqs. (fl"6|) and IT7|) . At and u will depend upon the time according to 

ATW = A T (0) v /l + 10^/A|(0)i, u{t) = ^=^B== (20) 

V y/1 + 10l//A r (0) t 

The corresponding values of a\ and ct2 are determined substituting these latter expressions into Eqs. 
(fl"8f and (fl9|) . respectively, i.e. 



_ A T d ( tu \ _ 1 
01 = Vdt I At" J ~ l + \Qvl\ 2 {0)t 

(21) 

_ 1 d\ T _ 5 
2 u dt R 

The coefficient ai decreases with the time and for t — > 00, ai — * 0, whereas 02 remains constant. 
Therefore, for £ — > 00, one obtains the self-similar correlation function /(f), which does not depend on 
the initial condition and that obeys to the following non-linear ordinary differential equation 



5 df 1-f df 2 (d?f idf\ 10 (22) 

R dr \ 2 dr R \ dr 2 r dr J R 



The first term of Eq. (|22p represents the time derivative of / and does not influence the mechanism 
of energy cascade. In line with von Karman j2], [4], we search the self-similar solutions over the whole 
range of f, with the exception of the dimensionless distances whose order magnitude exceed R. This 
corresponds to assume the self-similarity for all the frequencies of the energy spectrum, but for the 
lowest ones [2], 0]. Accordingly, the self-similar solutions of the von Karman-Howarth equation are 
also steady solutions, thus the first term of Eq. (|22|) can be neglected with respect to the other ones, 
i.e. 



The boundary conditions of Eq. (|23p are from Eqs. J2J) [2 

/(0) = 1,^ = (24) 
dr 

lim /(f) = (25) 

r — >oo 

Since the solutions / 6 C 2 (0, 00) exponentially tend to zero when r — * 00, and At is considered to be 
an assigned quantity, the boundary condition ([25]) can be replaced by the following condition in the 
origin 
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Therefore, the boundary problem represented by Eqs. (|23|) . (|24|) and (|25l) . corresponds to the 
following initial condition problem written in the normal form 

df 

(27) 




the initial condition of which is 

/(0) = 1,F(0)=0 (28) 

This ordinary differential system arises from Eq. (|23p . where according to Eq. (|26p . one must take into 
account that 

Hm m = 

?— >o f f^o df 



l im ^ = li m ^P=-l (29) 



4 Results and Discussion 

In this section, the solutions of the system (|27[) are qualitatively studied and numerically calculated 
for several values of the Reynolds numbers, and the statistics of the velocity difference is investigated 
through the analysis seen in the section [5] (Eqs. ([?))- (fTB")) ). 



The analysis of Eqs. (|27|) shows that, for large f, / exponentially decreases, thus all the integral 
scales of / are finite quantities and the energy spectrum is a definite quantity whose integral over the 
Fourier space gives the turbulent kinetic energy. 

It is worth to remark that, where K(r) is about constant, / behaves like 

/-Is f 2/3 (30) 

as a result, the energy spectrum varies according to the Kolmogorov law k -5 / 3 and the interval where 
this happens, identifies the inertial subrange of Kolmogorov. 

Furthermore, the solutions / G C 6 (0,oo), of the system (|27 |) -([28 |) . in the vicinity of the origin read 

as 

/(r) = l-y + ^^r 4 + 0(f 6 ) (31) 

The more rapid variations of / and of its derivative occur in close proximity of the origin 3 and these 
are responsible for the behavior of E(k) at high wavenumbers. Thus, we assume that the minimum 
length scale of / is due to the second and the third term of Eq. (f3Tj) . Under this scale the effects of the 
inertia and pressure forces are negligible with respect to the viscous forces. This scale is determined as 
two times the separation distance f > where the first derivative of — f 2 /2 + (10 + i?)/112 f 4 vanishes. 
Then, the corresponding wave-number (minimum scale) is related to R through the relationship 



_ 2ir 1Q + R 

For higher wavenumbers, in the dissipative range, the viscosity influences / in such a way that E(k) 
decreases much more rapidly with respect to the inertial subrange. In the dissipative interval, E{n) 
roughly coincides with the Fourier transformation of / ~ 1 — f 2 /2 ~ 1/(1 + ^ 2 /2), and, for this reason, 
one could expect that the energy spectrum behaves like [5] 

E(k) w exp(-aK) (33) 

for k > K maxi where a = 0(1) > is a proper parameter which depends upon the Reynolds num- 
ber. The value of n max , calculated following Eq. (|32p, should indicate the order of magnitude of the 
separation wavenumber between the Kolmogorov inertial subrange and the dissipative interval. 



Now, several numerical solutions of Eqs. (|2T[) were calculated for different values of the Taylor-scale 
Reynolds numbers by means of the fourth-order Runge-Kutta scheme of integration. 
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Fig. 1 Longitudinal correlation function for different Taylor-Scale Reynolds numbers. 
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Fig. 2 Longitudinal triple correlation function for different Taylor-Scale Reynolds numbers. 

The cases here analyzed correspond to R = 100, 200, 300, 400, 500 and 600. The fixed step size of 
the integrator scheme is selected on the basis of the asymptotic stability condition of Eq. (|27|) . This is 
Af = V2/R [6] and provides a fairly accurate description of the correlation function at small scales. 

Figures [T] and [2] show the numerical solutions of Eqs. ([27]) . where double and triple longitudinal 
correlation functions are both in terms of f, for the different values of R. Due to the mechanism 
of energy cascade, represented by the expression of K(r) (see Eq. the tail of / rises with R 

accordingly to Eq. (|3"Tj) and this means that, for an assigned value of At, all the integral scales of / rise 
with the Reynolds number. The triple correlation function is shown in terms of f in Fig. [2] Because 
of energy cascade, and according to Eq. k decays more slowly than /, for f — > oo. It is apparent 
that assigned variations of k correspond to variations of f whose size increases with R. The maximum 
of \k\ gives the entity of the mechanism of energy cascade. This is slightly less than 0.05 and agrees 
quite well with the numerous data of the literature (see [3] and Refs. therein) . 

Figures [3] and S] show the plots of E(k) and T(k) calculated with Eq. (fl"4"|) , for the same Reynolds 
numbers. As the consequence of the mathematical properties of /, the energy spectrum behaves like 
E(k) — 0(n A ) in proximity of the origin, and after a maximum is about parallel to the —5/3 Kol- 
mogorov law (dashed line in Fig. [3]) in a given interval of the wave-numbers. This interval defines the 
inertial range of Kolmogorov, and its size increases with R. This arises from the fact that there exists a 
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Fig. 3 Turbulent Energy Spectrum for different Taylor-Scale Reynolds numbers. 
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Fig. 4 "Transfer function T(re)" for several Taylor-Scale Reynolds numbers. 

spatial interval where K(r) is about constant and, as seen, this determines that /—1m r 2 / 3 . At higher 
wave-numbers, the energy spectrum decreases more rapidly with respect to the Kolmogorov subrange 
and its variations almost agree with Eq. ((33)) . The constant a of Eq. (|33[) varies from about 0.53 to 0.26 
when R assumes the values of 100 and 600, respectively. The wavenumber of separation between the 
two regions changes with the Reynolds number and follows the variations previously determined with 
Eq. (|32|) . Since K does not modify the kinetic energy of the flow, according to Eq. ([4}, the integral of 
T(re) over re results to be identically equal to zero at all the Reynolds numbers. 

In the Figs. [5] and [6l skewness and flatness of Au r are shown in terms of r for the same values of the 
Reynolds numbers. The skewness H 3 is first calculated according to Eq. §5$ and thereafter the flatness 
_ff 4 has been determined using Eq. ([5]). For a given value of R, \H 3 \ starts from 3/7 at the origin, then 
decreases to small values, while H4 starts from values quite greater than 3 at r = 0, then reaches the 
value of 3. In line with Eq. ((4]), since / decays more rapidly than k as r — > 00, H4 goes to 3 more 
fastly than H3 tending to zero. Although #3(0) does not depend upon R, H^(f) is a rising function of 
the Reynolds number and, in any case, the intermittency of Au r increases with R according to Eqs. 

©-CP- 
Next, the Kolmogorov function Q(r) and the Kolmogorov constant C, are determined using the 

previous results. According to the theory, the Kolmogorov function, defined as 

= jA^n (34) 

re 
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Fig. 5 Skewness of Au r at different Taylor-Scale Reynolds numbers. 
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Fig. 6 Flatness of Au r at different Taylor-Scale Reynolds numbers. 



is constant with respect to r, and is equal to 4/5 as long as r/X T = 0(1), where e = —3/2 du 2 /dt 
is the rate of the energy of dissipation. In Fig. [3 Q(r) is calculated through the skewness of velocity 
difference and is shown in terms of f. This exhibits a maximum for f = 0(1) and quite small variations 
for higher f , as the Reynolds number increases. These behavior is the consequence of the aforementioned 
variations of H^{r) with R and r. The maximum of Q rises with R and seems to tend toward the limit 
4/5 prescribed by the Kolmogorov theory. In Fig. [8l the maximum of the Kolmogorov function is in 
term of R. These data ("x" symbols) are shown in comparison with those of Ref. [7] (circular filled 
symbols and continuous line). There, the Kolmogorov equation with a forcing term, is compared with 
experimental measurements, for helium gas at low temperature. It is apparent that, for assigned _R, 
the corresponding values of Q max , calculated with Eq. (f34|) , are less than those of Ref. [7], and such 
difference varies with R with an average percentage value of about 12%. This disagreement could be 
due to the self-similarity here assumed or to possible differences in the estimation of the Taylor-scale. 
Specifically, At is here analytically calculated or assumed, whereas in Ref. [7], it is determined through 
measurement of dissipation rate. Nevertheless, the two set of data can be considered comparable, since 
the two variations almost exhibit the same trend. 

For what concerns the Kolmogorov constant C, it is defined by E(n) ps Ce 2 ^ / k 5 ^ 3 , and is here 
calculated as 



E(k)k 5 / 3 
G = max -7- — 

kG(0,oo) e 1 !* 



(35) 
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Fig. 7 Kolmogorov function for several Taylor-Scale Reynolds numbers. 
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Fig. 8 Maximum of the Kolmogorov Function in terms of the Taylor-Scale Reynolds number. These data are 
from Ref. [7]. The symbols "X" are for the present result. 



Table 1 Kolmogorov constant for different Taylor-Scale Reynolds number. 



R 


c 


100 


1.8860 


200 


1.9451 


300 


1.9704 


400 


1.9847 


500 


1.9940 


600 


2.0005 



In the tableQ] the Kolmogorov constant is shown in terms of the same Reynolds numbers. The obtained 
values of C increase with the Reynolds number and are in good agreement with the numerical and 
experimental values known from the various literature [8] , [9] , |10j . 

The spatial structure of Au r , expressed by Eq. ([7|), is also studied with the previous results. Ac- 
cording to the various works [llj . [12] . [13] , Au r behaves quite similarly to a multifractal system, where 
Au r obeys to a law of the kind Au r (r) w r q in which q is a fluctuating exponent. This implies that 
the statistical moments of Au r {r) are expressed through different scaling exponents C( n ) whose values 
depend on the moment order n, i.e. 



((Au r ) n ) = A n r^ n > (36) 

In order to calculate these exponents, the statistical moments of Au r are first calculated using Eqs. 
([8]) for several separation distances. Figure [9] shows the evolution of the statistical moments of Au r 
in terms of f , in the case of R = 600. The scaling exponents of Eq. (f3"6"|) are identified through a best 
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Fig. 9 Statistical moments of Au r in terms of the separation distance, for i?=600. 
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Fig. 10 Scaling expon ent s of Au r for several i?. Solid symbols are for the present data. Dashed line is for 
Kolm ogo rov K41 data [14]. Dotted line is for Kolmogorov K62 data [11]. Continuous line is for She-Leveque 
data [13] 



fitting procedure, in the intervals {f 1^2)1 where the endpoints r\ and fi have to be determined. The 
calculation of C( n ) an d A n is carried out through a minimum square method which, for each moment 
order, is applied to the following optimization problem 



J„(C(n),A„) 



\((Au r ) n ) - A n r^) 2 dr 



1,2, 



(37) 



where (((Au r ) n )) are calculated with Eqs. ([8]), fi is assumed to be equal to 0.1, whereas is taken 
in such a way that £(3) = 1 for all the Reynolds numbers. The so obtained scaling exponents are 
shown in Table ^ in terms of the Taylor scale Reynolds number, whereas in Fig. [TU] (solid symbols) 
these exponents are compared with those of the Kolmogorov theories K41 [14] (dashed line) and K62 
[TT] (dotted line), and with the scaling exponents calculated by She-Leveque [13] (continuous curve). 
Near the origin £(n) ~ n/3, and in general the values of £(n) are in good agreement with the She- 
Leveque results. In particular the scaling exponents here calculated are lightly greater than those by 
She-Leveque for n > 8. 

According to the present analysis, these peculiar laws C — C( n ) which make Au r similar to a 
multifractal system, are the consequence of the combined effect of the quadratic terms into Eq. (J7J 
and of the functions K and k calculated through Eq. (12"?) . 
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Table 2 Scaling exponents of the longitudinal velocity difference for several Taylor-Scale Reynolds number. 



R 


100 


200 


300 


400 


500 


600 


C(l) 


0.35 


0.35 


0.35 


0.35 


0.35 


0.35 


C(2) 


0.70 


0.71 


0.71 


0.71 


0.71 


0.71 


C(3) 


1.00 


1.00 


1.00 


1.00 


1.00 


1.00 


C(4) 


1.30 


1.29 


1.29 


1.29 


1.29 


1.29 


C(5) 


1.56 


1.55 


1.55 


1.55 


1.55 


1.55 


C(6) 


1.82 


1.81 


1.81 


1.81 


1.81 


1.81 


C(7) 


2.06 


2.05 


2.05 


2.04 


2.04 


2.05 


C(8 


2.31 


2.28 


2.28 


2.28 


2.28 


2.28 


C(9) 


2.53 


2.50 


2.50 


2.50 


2.50 


2.51 


C(io) 


2.76 


2.72 


2.73 


2.72 


2.72 


2.73 


C(n) 


2.97 


2.93 


2.93 


2.93 


2.93 


2.94 


C(12) 


3.18 


3.14 


3.14 


3.13 


3.13 


3.15 


C(13) 


3.39 


3.33 


3.34 


3.33 


3.33 


3.35 


C(14) 


3.59 


3.53 


3.54 


3.53 


3.53 


3.55 


C(15) 


3.79 


3.73 


3.73 


3.72 


3.73 


3.75 



The PDFs of du r /df can be formally determined with Eqs. (fl~3]) and ([7]). Specifically, these PDF 
are calculated with several direct simulations, where the sequences of the variables £, rj and £ are first 
determined by a gaussian random numbers generator. The distribution function is then calculated 
through the statistical elaboration of the data obtained with Eq. ((7]). The results are shown in Fig. 
[TTa and [TTb in terms of the dimcnsionless abscissa 

du r /df 



((du r /dr)y/ 2 

These distribution functions are normalized, in order that their standard deviations are equal to the 
unity. The figure represents the PDF for the several R, and the dashed curve represents the gaussian 
distribution functions. In particular, Fig. [TTb shows the enlarged region of Fig. [TTb . where 5 < s < 8. 
According to Eq. ([7|) , although the skewness of du r /df does not depend on R, the tails of PDFs change 
with R in such a way that the intermittency of du r /df rises with the Reynolds number. This increasing 
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Fig. 11 Probability distribution functions of the longitudinal velocity derivative for the different Taylor-Scale 
Reynolds numbers 
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intermittency, caused by the quadratic terms appearing into Eq. ([7|) , is the result of the constancy of 
the skewness of du r /df. 

5 Conclusions 

The obtained self-similar solutions of the von Karman-Howarth equation with the proposed closure, 
and the corresponding characteristics of the fully developed turbulence are shown to be in very good 
agreement with the various properties of the turbulent flow from several points of view. 
In particular: 

— The energy spectrum follows the Kolmogorov law in a range of wave-numbers whose size increases 
with the Reynolds number, whereas for higher wave-numbers, it diminishes according to an expo- 
nential law. 

— As the consequence of the skewness of velocity difference, the Kolmogorov function exhibits a 
maximum and relatively small variations in proximity of r = O(Xt). This maximum value rises 
with the Reynolds number and seems to tend toward the limit 4/5, prescribed by the Kolmogorov 
theory. 

— The Kolmogorov constant moderately rises with the Reynolds number with an average value around 
to 1.95 when R varies from 100 to 600. 

— The scaling exponents of the moments of velocity difference are calculated through a best fitting 
procedure in an opportune range of the separation distance. The values of these exponents are in 
good agreement with the results known from the literature. 

These results represent a further test of the analysis presented in Ref. [Ij which adequately describes 
many of the properties of the isotropic turbulence. 
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